function hfile2 = glm_PSCMaps(hfile)
% GLM::PSCMaps  - calculate PSC maps for single study GLM
%
% FORMAT:       pscmap = glm.PSCMaps;
%
% Output fields:
%
%       pscmap      VMP object

% Version:  v0.7b
% Build:    7083022
% Date:     Aug-30 2007, 10:39 PM CEST
% Author:   Jochen Weber, Brain Innovation, B.V., Maastricht, NL
% URL/Info: http://wiki.brainvoyager.com/BVQXtools

% argument check
if nargin < 1 || ...
    numel(hfile) ~= 1 || ...
   ~isBVQXfile(hfile, 'glm')
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid call to %s.', ...
        mfilename ...
    );
end
bc = bvqxfile_getcont(hfile.L);
if bc.NrOfStudies ~= 1 || ...
    bc.ProjectType ~= 1 || ...
    bc.TransformationType == 1
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid call to %s.', ...
        mfilename ...
    );
end
bbox = aft_BoundingBox(hfile);

% create VMP
hfile2 = newnatresvmp(bbox.BBox, bc.Resolution, 11 * ones(1, bc.NrOfPredictors - 1));
bc2 = bvqxfile_getcont(hfile2.L);

% calculate maps
for pc = 1:(bc.NrOfPredictors - 1)
    bc2.Map(pc).LowerThreshold = 0.5;
    bc2.Map(pc).UpperThreshold = 5;
    bc2.Map(pc).DF1 = bc.NrOfTimePoints - bc.NrOfPredictors;
    bc2.Map(pc).DF2 = 0;
    bc2.Map(pc).BonferroniValue = bc.NrOfVoxelsForBonfCorrection;
    bc2.Map(pc).Name = sprintf('PSC map for %s', bc.Predictor(pc).Name2);
    mapv = 100 * bc.GLMData.BetaMaps(:, :, :, pc) ./ ...
        bc.GLMData.BetaMaps(:, :, :, end);
    mapv(isinf(mapv) | isnan(mapv)) = 0;
    bc2.Map(pc).VMPData = mapv;
end

% set back
bvqxfile_setcont(hfile2.L, bc2);
